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ABSTRACT 

We present a model in which planetesimal disks are built from the combination of planetesimal 
formation and accretion of radially drifting pebbles onto existing planetesimals. In this model, the rate 
of accretion of pebbles onto planetesimals quickly outpaces the rate of direct planetesimal formation in 
the inner disk. This allows for the formation of a high mass inner disk without the need for enhanced 
planetesimal formation or a massive protoplanetary disk. Our proposed mechanism for planetesimal 
disk growth does not require any special conditions to operate. Consequently, we expect that high 
mass planetesimal disks form naturally in nearly all systems. The extent of this growth is controlled 
by the total mass in pebbles that drifts through the inner disk. Anything that reduces the rate 
or duration of pebble delivery will correspondingly reduce the final mass of the planetesimal disk. 
Therefore, we expect that low mass stars (with less massive protoplanetary disks), low metallicity 
stars and stars with giant planets should all grow less massive planetesimal disks. The evolution 
of planetesimal disks into planetary systems remains a mystery. However, we argue that late stage 
planet formation models should begin with a massive disk. This reinforces the idea that massive and 
compact planetary systems could form in situ but does not exclude the possibility that significant 
migration occurs post-planet formation. 

Subject headings: planets and satellites: formation, protoplanetary disks 


1. INTRODUCTION 

One of the most interesting findings of the Kepler mis¬ 
sion is the prevalence of super-Earth and Neptune sized 
planets on short period orbits. Many of these planets 
reside in compact multi-planet systems often referred to 
as STIPs (Systems with Tightly packed Inner Planets). 
Although such planets are lacking in the solar system, 
their ubiquity in the Kepler data indicates they are a 
n atural outcome of th e planet formation process. 

iChiang fc Laughllnl (.2013:) found that the minimum 
surface density of material needed to form many of the 
Kepler systems is at least 5-10 times higher than the sur¬ 
face density of the minimum mass solar nebula (MMSN) 
extrapolated inward of Mercury’s orbit. This high con¬ 
centration of mass within the inner 1 AU or so of these 
systems is a strong indication that the transport of solid 
material in protoplanetary disks is an important aspect 
of planet formation. The difficulty lies in determining 
how and when this material is transported. 

One hypothesis is that planets form first and then mi¬ 
grate to their observed locations. In this scenario plan¬ 
ets would form farther out in the disk and migrate in¬ 
wards due to gravitational torque s exerted on them by 
the disk (see Kiev fe Nelsonl 120121 for review). N-body 
simulations of planet formation that include migration 
forces have successfully produced systems of, ty pically 2- 
5, planets on close in or bits (e.g lTerauem fc Papaloizoul 
120071 : lHands et al.ll2014l f. However, the formation with 
migration scenario is not without its problems. Most no¬ 
tably, these simulations produce an over-abundance of 
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neighboring planets that are at or near mean motion res¬ 
onances when compared to the Kepler sample. It is pos¬ 
sible that additional forces that act after the planets have 
reached thei r final position may remove planets from res¬ 
onan ce (e.g. iLithwick fc Wull2012l : iBatvgin fc Morbidcllil 
[ 20131) . Another prediction of the migration scenario is 
that planets form beyond the ice line and consequently 
will be very water rich. Observations of planet masses 
and radii find that of the planets with masses known to 
better than 20%, the smaller planets (those with radii 
less that 1.7 R®) have den sities that are consist ent with 
an Earth-like composition dDressing ctalj |2015|L _ 

Ano ther possibility, promoted bv lHansen fc Murravl 
(I 20 H , is that solid material migrated through the disk 
before planet formation and formed a massive planetesi¬ 
mal disk. Planets then formed from this disk in a similar 
manner as we believe the terrestr i al pla nets of the so¬ 
lar system did. lHansen fc Murra^ (|2012H simulated the 
in situ formation of close in super-Earths and Neptunes 
from a massive disk and found good agreement between 
the observed distribution of planets and their simulated 
systems both in terms of the mass-period distribution 
and more detailed statistics. The major downside to this 
hypothesis is the omission of the physics of migration, 
which is believed to be important for bodies larger than 
Earth. 

One of the difficulties in determining which, if either, 
of these scenarios is correct is our lack of knowledge of 
the initial conditions of planet formation. Specifically, 
the surface density of solids throughout the disk directly 
preceding planet formation is poorly known. This uncer¬ 
tainty adds a huge amount of freedom to planet forma- 
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tion models which is essentially unbounded. 

Our poor knowledge of the surface density distribution 
in planetesimal disks is understandable considering the 
lack of observational constraints. What information we 
do have comes from the remnants of the solar system’s 
planetesimal disk. However, we can still make predic¬ 
tions of what we expect planetesimal disks to look like 
which can offer some constraint to planet formation mod¬ 
els. 

One crucial piece of information gleaned from mete¬ 
orites is that planetesimal disks are built u p over time - 
over the course of about 2-3 million years ([Amelin et al.l 
I2002D . With this knowledge, we can think of planetesimal 
disk formation as an evolutionary sequence and not just 
some scaling of the gas/dust content of a protoplanetary 
disk at some point in time. In this context, the growth of 
a planetesimal disk depends on the rates of planetesimal 
formation, growth, destruction and migration. Integrat¬ 
ing these rates over time, throughout the disk, provides 
us with the surface density of planetesimals to use as 
initial conditions for planet formation models. 

In order to understand the growth of the planetesi¬ 
mal disk, we need to know what the reservoir of mate¬ 
rial is that is fueling this growth. Models of dust evo¬ 
lution in protoplanetary disks predict rapid growth of 
grains to millimeter and centimeter-size pebbles, (e.g. 
iBrauer et al.l 120081 iBirnstiel et al.1 l2010l) . This predic- 
tion is supported by millimeter observations of disks 
(jTesti et al.ll20l4 and references therein). At this size, 
pebbles should experience fast radial drift due to the 
headwind they feel orbiting in the sub-Keplerian gas disk. 
This constant supply of drifting pebbles is a likely reser¬ 
voir for planetesimal growth. 

Assuming that radially drifting pebbles do indeed sup¬ 
ply the inner disk with solids, the growth of the planetes¬ 
imal disk boils down to how a nd where this mat e rial is 
halted from further migration. iChatteriee fe Tanl (|2014f ) 
proposed that migration of pebbles will be halted by a 
pressure bump in the inner disk associated with the inner 
edge of an MRI “dead-zone”. This would cause pebbles 
to collect in a dense ring at which point they would either 
become gravitationally unstable and collapse into proto¬ 
planets or they would undergo planet formation via core 
accretion. As this dead-zone boundary retreats, away 
from the star subsequent planet formation could occur, 
thus building STIPs from the inside out. An alterna¬ 
tive solution is that pebbles simply grow fast enough and 
large enough in the in ner disk that further radial drift is 
halted. iBolev ct~ahl (12014) suggest that this rapid growth 
may be plausible due to a high concentration of solids 
in the inner disk and very efficient sticking of partially 
molten pebbles. 

We propose an alternative means to stop pebbles from 
further migration. We require that some seed planetesi¬ 
mals form in the inner disk. As we will show, the initial 
seed population can be very small. The rate of pebble 
accretion onto these existing planetesimals quickly sur¬ 
passes the rate of new planetesimal formation in the inner 
disk such that the majority of mass in the final plan¬ 
etesimal disk is from the accretion of pebbles, not from 
the direct formation of new planetesimals. The appeal 
of this model is the generality of our assumptions. The 
formation of planetesimals almost certainly occurs as evi¬ 
denced by the small bodies of the solar system which are 


presumably remnants of a planetesimal disk. It seems 
reasonable to extend this assumption to the inner disk 
as most models of planetesimal formation all ow for the 
form a tion of planetesimals in this region (e.g . I Chambers! 
1201(1 : IBolev et~ahl 12014 : iCarrera et al.l l2015[ ). However, 
uncertainties in the planetesimal formation process limit 
our ability to evaluate the validity of this assumption. 
The assumption of a supply of pebbles from the outer 
disk is consistent with both models and observation. Ad¬ 
ditionally, our model does not require the presence of a 
pressure bump or near perfect sticking between pebbles. 
Instead, high mass inner disks are a natural consequence 
of global accretion physics. 

2. MODEL 

The growth of a planetesimal can be described as: 

Mp —7 lb Ppeb^reh (1) 

where b is the maximum distance a planetesimal can ac¬ 
crete a pebble from, p peb is the volume density of peb¬ 
bles near the planetesimal and v re i is the relative velocity 
between the pebbles and the planetesimal. In the case 
where the accretion radius, 6, is larger than the scale 
height of the pebble disk, this can be written as: 

M p = 2bY, peb v re i, (2) 

where T, peb is the surface density of pebbles in the disk. 
In the following sections we describe how the density of 
pebbles in the disk and the accretion rate of these pebbles 
are determined. 

2.1. Protoplanetary disk profile 

The gas density throughout the disk is a key piece of in¬ 
formation that is necessary for other parts of our model. 
We use a simple prescription for the gas surface density 
profile: 

^.=0(f u y r . o) 

We choose p = 1, whi ch is predicted for a steady- state 
viscous accretion disk (iLvnden-Bcll &; Pringldll974 ) and 
is also consistent with the range o f values fitted to sub- 
millimeter obser vations of disks (iA ndrcws~ fc Williams! 
I2007T) . Following lLambrcchts fc Johansenl (120141 . we set 
/3 = /3oexp(—f/3Myr) to account for the dissipation of 
the disk over time. The gas volume density is related to 
the surface density by p gas ~ T, gas /2H, where the scale 
height is given by: 

H / r N 1 / 4 

7“' 0033 Gw) ' <4) 

2.2. Pebble surface density 

In our model, we assume that pebbles are continuously 
supplied from the outer disk and that there is no contri¬ 
bution to the pebble flux from within the simulation re¬ 
gion (< 5AU). The only sink for pebbles is by accretion 
onto planetesimals. The mass flux of pebbles through 
the disk at a given location is then: 

Fpeb(r') = F sup ) Mp = 27rrv r Yjp eb: ( 5 ) 

r p >r 

where F sup is the supply rate of pebbles to the simula¬ 
tion region and the sum is over the accretion rates of all 
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planetesimals outside of r. Inverting this equation gives 
us the local surface density of pebbles, 


frequency and the gas drag timescale, Td, is 


Zipeb — 


F, 


peb 


2nrv r 


( 6 ) 


The timescale for pebbles to drift into the star from the 
outer disk tends to be much shorter than the timescale 
for pebble growth. Consequently, pebble growth can be 
considered the limiting factor in setting F sup and Equa¬ 
tion [Hl_cm3j3e_considered_th^nstajitaneous surface den¬ 
sity. ^a^brechti^Johanie^ (120141) have developed a 
model to describe the growth of pebbles and the evolu¬ 
tion of the inward pebble mass flux over time, which we 
adopt for our work. For brevity, we do not repeat their 
derivation but instead, give their result. Assuming a gas 
disk profile as shown in Eq. [3] they found the supply of 
pebbles over time to be 


F = 9 5 x 10- 5 (—?—\ ( Ml )'' 3 ( FlY ' 3 

sup ' \500gcm- 2 ) \M Q ) ^ 0 . 01 / 


^500 gem 2 J \Mq J 
f t \~ 1/3 

,IoV j Meyr ~' 


(7) 


where Zq is the dust to gas mass ratio and /3 the gas 
surface density at 1AU, decays exponentially with time 
as in section 2.1. 

At some point the supply of pebbles will cease. 
This will occur either when the disk has dissipated 
or when the pebble productio n front, as described in 
lLambrechts fc Johansen! (j2014l ), moves beyond the outer 
radius of the protoplanetary disk. Assu min g the typica l 
exten t of a disk is around 200 AU ((Andrews fc Williams! 
120071b the pebble supply would cease after 1-2 million 
years. We choose our simulation duration to be one 
million years as it is shorter than both the dissipation 
timescale and the pebble depletion timescale. This en¬ 
sures that we do not supply more pebble mass to the 
inner disk than is available in the whole protoplanetary 
disk. It is possible that pebble accretion will occur longer 
than our simulation duration, but as we show later, this 
is similar to increasing the pebble production rate, which 
we do test. 

The radial drift velocity of the pebbles is given by 


v r = — 2- 


St 


St 2 + 1 


V^kep 


( 8 ) 


ljWeidenschillingtl977tlNakagawa et al.lll986T) . where v kep 
is the orbital Keplerian velocity, rj describes the pressure 
support on a body and can be estimated as 


rj = 0.0018 


Td = 


PmS 

Pgas 

4 p m S 2 


QPgasCs^ 1 
6 p m S 
PgasP^kep 


9A 

s < — (Epstein drag) 


9A 27c s A , . 

— < s < — - (Stokes drag) (10) 

4 2r)v kep 


27c s A 

s > - (Quadratic drag). 

2 pVkep 


p m is the material density of the pebble, s is the size of 
the pebble, c s is the local so und speed of the gas a nd 
A = 2 x 10 ~ 9 /p gas (cgs units; iNakagawa et al.l 1198(f ) is 
the mean free path of the gas. 

Substitution of Eqs. [7] and [8] into Eq. [G] yields the sur¬ 
face density of pebbles throughout the disk over time. 
The pebble volume density can be approximated by di¬ 
viding the surface density by 2 times the scale height of 
the pebbles: 


Hpeb — I! 


a 


■St 


1/2 


( 11 ) 


(jYoudin fe Lithwicki 12007T) , where a is the turbu- 
lent diffusion parameter for a viscous accretion disk 
dShakura fc Sunvacvlll973lh 

2.3. Pebble Accretion Rate 

To obtain the accretion rates of planetesimals, we need 
to know their collisional cross-section. In the simplest 
case, this is the geometrical cross-section. In reality, 
other factors act to enhance the cross-section above the 
geometrical limit such as gravitational focusing and gas 
drag. Gas drag, in particular, is critically important 
when considering the accretion of small pebbles because 
its affect on the trajectories of pebbles is of compara¬ 
ble m agnitude as the affect of gravity. lOrmel fc Klahrl 
• 201( 1: calculated the trajectories of small bodies near 
planetesimals and found that the combined influence 
of gravity and gas drag on the small bodies can in¬ 
crease a planetesimal’s accretional cross-section signif¬ 
icantly - up to a large fraction of the planetesimal’s 
hill radius. Subsequent studies have employed this re¬ 
sult to help explain the short timescales that are neces¬ 
sary to f orm gas giant cores before the dispersal o f the 
gas disk (lLambrechts fc Joh ansen 20 liS lChamberdl20l4 


lLambrechts fc Johansenll2014t iKretke fc Levison 2014). 

In our work, we apply the results of lOrmel fc Klahrl 

(l201(tl to the accretion of pebbles in the inner disk. We 
direct the reader to their paper for details of these calcu¬ 
lations. In particular, their Table 2 provides a summary 
of the steps we took to determine the collisional cross- 
section of planetesimals. 

The remaining term that is needed to determine the 
accretion rate of pebbles is the relative velocity between 
the planetesimal and the pebbles. We take this to be the 
larger of either the radial drift velocity (Eq. [HJ or the 
velocity imparted upon the pebbles due to the turbulent 
motions of the gas, 


(12) 


/ r \ !/ 2 



(w) ■ 

(9) 

v'turh - ac 2 s x mm ( 2St, x + gt j 


St = Tdll is the Stokes number where Q is the orbital (lOrmel fc Cuzzil[2007l : lChambersll2014f) . 





















































4 


TABLE 1 

Simulation paramters 


Parameter 

Default value 

Description 

Ot-v 

0.001 

Viscosity parameter 

s 

1 cm 

Pebble size 

Ppeb 

2 gem -3 

Pebble material density 

R P 

100-500 km 

Planetesimal size 

Pp 

3 gem -3 

Planetesimal material density 


500 gcm~ 2 

Initial gas surface density at 1 AU 

So 

1 gem -2 

Initial planetesimal surface density 
at 1 AU 

a 

-1.5 

Initial planetesimal surface density 
power law index 


The accretion rate is calculated first for the outermost 
planetesimal then proceeding inward to the innermost. 
The surface density of pebbles at a planetesimal location 
is calculated according to Equation^ The accretion rate 
is then calculated by inserting this value, the relative 
velocity between the planetesimal and pebbles and the 
calculated accretion radius into Equations [T] or [2] The 
inward pebble flux is then decreased by the accretion rate 
of that planetesimal in order to conserve mass and the 
process is repeated for the next planetesimal. 


2.4. Initial Conditions 

We begin the simulations with a low mass disk of plan- 
etesimals. Conceptually, the start of our simulations cor¬ 
responds to the point at which the rate of pebble accre¬ 
tion exceeds the rate of new planetesimal formation. Be¬ 
cause, the process of planetesimal formation is not well 
understood, the distribution of the first planetesimals is 
poorly constrained. For this reason, we test a range of 
initial planetesimal distributions with the following func¬ 
tion form: 



The locations of the planetesimals are determined by 
assuming a surface density profile for the planetesimal 
disk. Planetesimals, with some size range/distribution, 
are then placed in the disk so that they follow the as¬ 
sumed surface density profile. 

The expected sizes of newly born planetesimals de¬ 
pends on how they form. If they form via pairwise ac¬ 
cretion then a continuum of sizes from small to large 
would be expected. If instead they form from the grav¬ 
itational collapse of pebbles in either turb ulent concen¬ 
trations dCuzzi et al.ll2008HCham bcrs 2010|1 or streaming 
instabilities ( Johansen et al.ll2007fl then we would expect 
the initial size distribution of planetesimals to be larger. 
The size distribution of bodies in the asteroid belt also 
supports the idea that planetesi mals (at least in tha t re¬ 
gion of the disk) were born big (jBottke et al.l [2005h . In 
our simulations, we also assume that planetesimals form 
large. For most of our simulations, the initial planetes- 
imal sizes follow a distribution similar to the inferred 
initial size distribution of the asteroid belt dBottke et all 

f200l : 


N(> D) oc £> -4 ' 5 100 km < D < 500 km. (14) 


3. SIMULATIONS 
3.1. Baseline Simulation 
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Fig. 1.— Inward flux of pebbles throughout the disk at different 
times during the simulation. The reduction in pebble flux closer 
to the star is due to pebble accretion by planetesimals. 



Fig. 2.— Surface density of planetesimals throughout the disk 
at different times during the simulation. Red and blue lines corre¬ 
spond to the MMSN and 10 times the MMSN respectively. 

We begin by examining the general behavior of the 
accretion simulations for a baseline case. Our model pa¬ 
rameters and their values for this baseline case are listed 
in Table [TJ In the following sections we explore how the 
choice of parameter values affects the the growth of the 
planetesimal disk. 

A striking feature of this simulation is the high effi¬ 
ciency with which pebbles are accreted onto planetesi¬ 
mals. Figure [T] shows the inward flux of pebbles through 
the disk as a function of semi-major axis for a number of 
different times in the simulation. Except at the very be¬ 
ginning of the simulation, nearly all pebbles are accreted 
before they can make it through the disk. Wc will refer 
to the earliest time at which all pebbles are accreted as 
the time of total accretion. 

The pebble accretion rate is a strong function of loca¬ 
tion in the disk. The reason for this is twofold. The sur¬ 
face density of planetesimals increases closer to the star. 
Therefore, any particular pebble drifting in through the 
disk has a higher probability of being accreted by a plan¬ 
etesimal the closer it gets to the star. The second rea¬ 
son is that as pebbles drift inwards they are concentrated 
from a relatively large volume into one much smaller thus 
increasing the space density of pebbles dramatically at 
small radii. Planetesimals accreting from this concen¬ 
trated flux of pebbles will grow rapidly. Because the 
accretion rate is such a strong function of distance from 
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the star, the radial transition in the disk from minimal 
accretion to significant accretion is sharp. This edge of 
accretion moves out in the disk over time as planetesi- 
mals farther out in the disk grow large enough to begin 
accreting significantly from the pebble flux. Interior to 
this edge the accretion rate begins to drop off again. This 
decrease is due to the fact that most of the pebbles have 
already been accreted by planetesimals farther out and 
so the supply is much reduced. This forms an inner edge 
to the accretion. Accretion within this annulus accounts 
for the majority of pebble accretion throughout the disk 
at any given time. As the accretion annulus progresses 
outward over time, the planetesimal disk is built up from 
the inside out, leading to the evolution depicted in Figure 

m 

3.2. Initial Planetesimal Size 

The initial sizes of planetesimals are not known. Evi¬ 
dence from the asteroid belt suggests that, in that region 
of the disk, planetesimals were born large. For most of 
our simulations we begin with a distribution of planetes¬ 
imal sizes that are consistent wi th that inferred from the 
asteroid belt (jBottke et al.ll2005h . However, we do not 
know if this size distribution is accurate for the rest of 
the disk. Here we consider the effects of changing the 
initial sizes of the planetesimals. 

The importance of planetesimal size is best demon¬ 
strate d by examining the results of lOrmel fc Klahrl 
( 201( 1). In their work, they numerically integrated the 
trajectories of small particles near planetesimals and in¬ 
cluded the effects of gas drag. The combined effects of 
gravity and gas drag lead to a somewhat complicated 
dependence of the accretional cross-section on planetesi¬ 
mal size, pebble size and disk properties. In Figured! we 
calculate the pebble accretion rate, using the collisional 
cross-sections from (jOrmel fc Klahrl 12010 1, throughout 
the disk for disks composed of planetesimals of different 
sizes. The total pebble accretion rate in the inner disk 
depends weakly on the size of the planetesimals except 
at very large sizes at which point accretion becomes less 
efficient. Beyond ~ 1 AU, a disk composed of larger plan- 
etesimals accretes significantly faster than one composed 
of smaller planetesimals. 

With this insight, we can interpret the results of the 
accretion simulations with different initial planetesimal 
sizes. Figure Q] shows the final surface density of the 
planetesimal disk for simulations with four different ini¬ 
tial planetesimal sizes. As might be expected, disks with 
larger planetesimals accrete more pebbles farther out in 
the disk. This leads to a reduced pebble flux closer to 
the star, and, consequently, a less massive inner disk. 
Disks composed of larger planetesimals are also less ef¬ 
ficient at accreting planetesimals over the course of the 
full simulation. This is because planetesimals larger than 
~ 2000 km are less efficient (per unit mass) at accreting 
pebbles. The initially larger planetesimals will cross this 
size threshold sooner and therefore accrete less efficiently 
over their lifetime. 

3.3. Initial Planetesimal Distribution 

The initial placement and total mass of planetesimals 
in the disk are perhaps the most ill-constrained aspects 
of these simulations. In reality, the distribution of the 



Fig. 3.— Relative accretion rate of pebbles for disks composed 
of different sizes ( R p in units of km) of planetesimals. Each disk 
has the same total mass and radial surface density profile. 



Fig. 4.— Final surface density of the planetesimal disk in simula¬ 
tions starting with different initial planetesimal sizes ( R p in units 
of km). Red and blue lines correspond to the MMSN and 10 times 
the MMSN respectively. f a gives the fraction of pebbles that were 
accreted during the simulation. 

first planetesimals must depend on the details of plan¬ 
etesimal formation and how they change throughout the 
disk. Until we have a better understanding of planetes¬ 
imal formation we must base our initial conditions on 
other information. 

Observations of gas in protoplanetary disks provide in¬ 
formation on the distribution of mass p r ior to and dur ing 
planetesimal formation. lAndrcws fc Williamsl (12007D fit 
power laws to the gas surface density as a function of 
disk radius and found the power law index to range be¬ 
tween 0 and -2. On the other end, observations of plan¬ 
etary systems provide us with the final distribution of 
mass once the planet formation process has completed. 
On the whole, exoplanetary s ystems have a su rfa ce den - 
sity power law index of -1.5 (iChianer fc Laughhnll2013 i. 
consistent with the MMSN, but ind ividual systems have 
a wid e range of power law indices dRavmond fc Cossoul 
I2014T) . These constraints provide us with a starting 
point, but it is important to keep in mind that they are 
not strict constraints as the solid mass distribution in the 
intermediate stages of planet formation does not neces¬ 
sarily reflect the mass distribution in the initial and final 
stages. 

In most of our simulations, we began with the initial 
planetesimal distribution following a tenth of the MMSN. 
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Fig. 5. — Final surface density of the planetesimal disk in simu¬ 
lations starting with different initial total masses in planetesimals. 
So is the initial surface density of planetesimals at 1 AU in units 
of gcm~ 2 . Red and blue lines correspond to the MMSN and 10 
times the MMSN respectively. f a gives the fraction of pebbles that 
were accreted during the simulation. 

Here we remove this restriction and allow both the total 
mass in planetesimals and the slope of the mass distri¬ 
bution to vary. Figure [5] shows the final planetesimal 
surface density for three different simulations beginning 
with a planetesimal surface density distribution 0.001, 
0.01 and 0.1 times the MMSN. The effects of changing 
the initial amount of mass in the system are relatively 
limited. Simulations with initially lower mass planetesi¬ 
mal disks are less efficient at accreting pebbles over the 
course of the simulation. Nevertheless, even the very low 
mass disk accretes the majority of pebbles. This differ¬ 
ence can be attributed to the initially lower mass disk 
taking longer to reach the time of total accretion (time 
where all pebbles are accreted) and the relatively small 
magnitude of this difference can be attributed to the time 
of total accretion being much less than the length of the 
simulation. Similarly, changing the slope of the initial 
planetesimal surface density profile has little impact on 
the final surface density distribution of planetesimals (see 
Figure [6] 

An interesting feature is that the simulation with the 
initially lower mass disk ends up with a slightly higher 
surface density in the innermost regions than the other 
simulations. This results because the accretion annulus 
moves out in the disk more slowly resulting in a longer 
period of accretion in the inner disk. On the other hand, 
this results in a shorter period of accretion for the outer 
portions of the disk and a consequently lower final surface 
density in the outer disk. 

3.4. Pebble Size 

The size of the pebbles being accreted has a signifi¬ 
cant effect on the accretion rate of a planetesimal. The 
typical size of pebbles in the inner regions of disks is 
unknown, but we can test a range of values that are 
motivated from dust coagulation simulations. Such sim¬ 
ulations find t he dominant particle size to range from 
0.01-10 cm (e.g. lBirnstiel et aLil2010HBraucr etahll2008( l. 
Figure 0 demonstrates the difference in accretion rate 
throughout the disk if we vary the pebble size over this 
range. It is clear that the accretion rate’s dependance 
on pebble size is complicated and even more so if we 
consider that the size of the accreting planetesimal also 



Fig. 6.— Final surface density of the planetesimal disk in sim¬ 
ulations starting with initial planetesimal surface density profiles 
with different power law indices (a). Red and blue lines correspond 
to the MMSN and 10 times the MMSN respectively. f a gives the 
fraction of pebbles that were accreted during the simulation. 



Fig. 7.— Relative accretion rate of pebbles of different sizes (s 
in units of cm) throughout the disk by 200 km planetesimals. In 
this figure, unlike the simulations, all parts of the disk receive the 
same radial pebble flux. 

changes the pebble size dependance. 

When we vary the pebble size in our simulations, we 
see a corresponding change in the final planetesimal sur¬ 
face density (Figure [Sj. Somewhat counter-intuitively, 
the smaller pebbles, which tend to be more efficiently 
accreted, result in less massive inner disks. Once again 
it is easiest to think in terms of the movement of the ac¬ 
cretion annulus. In simulations with small pebbles, this 
annulus moves out in the disk rather rapidly so that the 
innermost regions do not have time to accrete a lot of 
mass. This is then compensated by more accretion far¬ 
ther out. In all cases the efficiency of pebble accretion for 
the disk as a whole is very efficient, it is just the distri¬ 
bution of this mass that changes as a function of pebble 
size. 

3.5. Turbulence 

Turbulence in the disk acts to increase the scale height 
of pebbles in the disk (Eq. fTTIl . An increase in scale 
height corresponds to a decrease in the volume density 
of pebbles and consequently a reduction in the efficiency 
of accretion. A reduction in accretion efficiency slows the 
outward movement of the accretion annulus resulting in 
more mass being accreted closer to the star (see Figure 
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Fig. 8. — Final surface density of the planetesimal disk in simu¬ 
lations with different pebble sizes (s in units of cm). Red and blue 
lines correspond to the MMSN and 10 times the MMSN respec¬ 
tively. f a gives the fraction of pebbles that were accreted during 
the simulation. 



Fig. 9. — Final surface density of the planetesimal disk in sim¬ 
ulations with different turbulent diffusion parameters (ol u ). Red 
and blue lines correspond to the MMSN and 10 times the MMSN 
respectively. f a gives the fraction of pebbles that were accreted 
during the simulation.. 

ED. Increasing the strength of turbulence does reduce the 
efficiency of pebble accretion, but not enough to prevent 
the formation of massive planetesimal disks. Instead it 
has more of an impact on where in the disk this mass is 
accreted. 

3.6. Gas Mass 

One of the more interesting parameters to examine is 
the total mass of the disk. All the previous parameters 
we examine, we varied to account for their uncertainty. 
In reality these parameters may be very similar in differ¬ 
ent disks. Disk mass, on the other hand, is observed to 
vary from disk to disk and consequently has some pre¬ 
dictive power with regards to the distribution of mass in 
planetary systems. 

Varying the mass of the disk in our model is effectively 
varying the supply rate of pebbles - a higher mass disk 
leads to a larger supply of pebbles. The dust to gas 
ratio also determines the pebble supply rate so these two 
parameters are mostly degenerate. We will restrict our 
discussion to the total disk mass, but the same effects 
can be achieved by varying the dust to gas ratio in a 
corresponding manner. 

The effect of changing the disk mass is relatively 



Fig. 10.— Final surface density of the planetesimal disk in sim¬ 
ulations with different initial gas disk masses. The parameter, /3o, 
is the initial gas surface density at 1 AU in units of g cm ~ 2 . Red 
and blue lines correspond to the MMSN and 10 times the MMSN 
respectively. f a gives the fraction of pebbles that were accreted 
during the simulation. 

straight-forward. A lower mass disk essentially slows 
down the whole accretion sequence. The accretion an¬ 
nulus moves out in the disk much more slowly. This 
results a lower mass planetesimal disk and a larger frac¬ 
tion of the accreted mass lying in the innermost regions 
of the disk (Figure [Till). Note the similarity between the 
evolutionary sequence of the growing planetesimal disk 
shown in Figure El and the different mass disks of Figure 
m In this sense, changing the disk mass has a similar 
effect to suddenly cutting off the supply of pebbles to the 
inner disk. 

3.7. Concurrent planetesimal formation and accretion 

A common feature in the results of many of the pre¬ 
vious simulations is a lack of mass in the planetesimal 
disk at more distant radii. This is a consequence of the 
longer planetesimal growth timescales and at these dis¬ 
tances. In all of these simulations, we began with some 
amount of mass in planetesimals and assumed that the 
pebble accretion rate was much larger than the rate of 
new planetesimal formation (and thus did not include 
any further planetesimal formation). However, this is 
clearly not the case for the outer disk where planetesi¬ 
mal formation must continue to be significant in order to 
produce enough mass to match, for example, the MMSN. 

To address this problem, we ran a simulation that con¬ 
tinues to form planetesimals throughout the course of the 
simulation. Because the rate of planetesimal formation 
in disks is unknown, we take the simplest case and as¬ 
sume that it is constant in time for the extent of the 
simulation. We start with an empty disk and at each 
timestep add planetesimals to the disk such that by the 
end of the simulation, were there no pebble accretion, 
the planetesimal disk would have a surface density pro¬ 
file corresponding to the MMSN. 

The results of this simulation can be seen in Figure 
□H Not surprisingly, the surface density profile at larger 
radii, where pebble accretion was negligible, matches the 
MMSN. In the inner disk, where pebble accretion was 
not negligible, the surface density is ten times higher 
than the MMSN, consistent with our previous simula¬ 
tions. This simulation is likely more realistic than the 
previous simulations because planetesimal formation and 









Fig. 11.— Final surface density of the planetesimal disk in a sim¬ 
ulation that included pebble accretion and planetesimal formation 
throughout the entire simulation. The model parameters are the 
same as those in Table 1 except that this simulation begins with 
no planetesimal disk. Red and blue lines correspond to the MMSN 
and 10 times the MMSN respectively. 

pebble accretion probably occur concurrently. However, 
we do not explore this scenario further because of the 
large uncertainties in the rate of planetesimal formation. 
The combination of planetesimal formation and pebble 
accretion throughout the lifetime of the disk is a realistic 
way to form a massive inner disk without causing an un¬ 
derdeveloped outer disk. This allows for the formation of 
STIPs and the formation of systems like the solar system 
from the same model, albeit with the requirement that 
some other m echanism eventually remove the i nne r disk 
material (e.g. iBatvgin fc Laughlid 120 1 5t IPu fc Wull2015l : 

lYolk fc Gladmanll2015f) 

4. DISCUSSION 

The most prominent result from our simulations is that 
the accretion of pebbles in the inner disk is an efficient 
process. Furthermore, this efficiency depends weakly on 
the model parameters. In all but one simulation (the 
very low mass protoplanetary disk simulation), the effi¬ 
ciency of pebble accretion was greater than 50% and in 
most cases much higher than that. Our results suggest 
that high mass and compact planetesimal disks are the 
expected outcome of planetesimal formation and growth 
in protoplanetary disks. 

4.1. Variation between planetary systems 

The processes described by our model are expected to 
operate in all protoplanetary disks. The observed distri¬ 
bution of planets shows a diversity of planetary system 
architectures. Our model naturally allows for this vari¬ 
ability of observed systems. 

The main source of differences between systems in our 
model is the extent of accretion. This is controlled by 
the duration of accretion and the rate of pebble delivery. 
Planetesimal disks that have more time to grow or receive 
a larger pebble flux are more massive and they are less 
centrally condensed than disks that have less time to 
grow or receive lower pebbles fluxes (see Figures [2] & 

m- 

Planetesimals will continue to accrete until the supply 
of pebbles is shut off. This may last as long as it takes 
for the disk to dissipate or in the event of giant planet 
formation, it may have a much shorter duration. A gi¬ 


ant planet in the outer disk will accrete the vast major¬ 
ity of pebbles drifting by thus significantly reducing the 
pebble flux to the inner disk. If giant planet formation 
occurs early on (within 1-2 million years), when pebble 
accretion is growing the inner disk, then we would ex¬ 
pect that, on average, systems containing giant planets 
would have a lower total mass in terrestrial planets and 
that the terrestrial planet mass would be concentrated 
in shorter period orbits. Individual systems containing 
giant planets would likely show a continuum of archi¬ 
tectures depending on the exact timing of giant planet 
formation. 

The total flux in pebbles to the inner disk depends 
on two main factors: the total mass of the disk and 
the solid fraction of the disk. _ Observations of disks 
(e.g. lAndrews fc Williamd I2007T) reveal a broad range 
of protoplanetary disk masses spanning several orders of 
magnitude. Similarly, t he span of metallicities f ound in 
planet hosting stars fe.g. lValenti fc Fischerll2005lf implies 
a broad range in dust to gas ratios in the disks from which 
they formed. These large ranges in disk mass and gas to 
dust ratios would lead to significant variations in the sup¬ 
ply rate of pebbles in different disks and consequently to 
variations in the structure of the planetesimal disk. 

Protoplanetary disk mass has been show n to corre¬ 
late with stellar type (jAndrews ct al.1 120131) . We might 
then expect planetesimal disks around larger stars to 
be more massive and less compact. However, the 
formation of giant planets complicates matters. Gi¬ 
ant planet formation has also bee n show n to corre¬ 
late with stellar mass ([Johnson et al.ll2007li . Therefore 
the expected trend with stellar mass may be washed 
out or removed due to the increased frequency of gi¬ 
ant planets. Similarly, giant planet occurrence rates 
are known to c o rrelate with host st a r metallicity (e.g. 
Gonzalej 11997!: [Santo s et all 1200 It iFischer & V alentil 

20051 iBuchhave et all 12014 iWang fc Fisched 1201511 no- 
tentially negating possible trends with metallicity. 

Our simulations also showed that changing other pa¬ 
rameters, such as the initial planetesimal sizes or the tur¬ 
bulent diffusion parameter, changes the resulting density 
distribution of the planetesimal disk. The range of pa¬ 
rameter values that we tested were meant to reflect their 
uncertainty and not their intrinsic variation between pro¬ 
toplanetary disks. However, some variation may still oc¬ 
cur between disks, especially in disks around different 
type stars leading to additional variability in planetesi¬ 
mal disk structure. 

In summary, the variation in protoplanetary disk 
masses and gas to dust ratios and the variation in the 
duration of pebble accretion give rise to variations in 
planetesimal disk structure. These factors affect the sur¬ 
face density structure of the planetesimal disk in similar 
ways and may, in some cases, counteract each other. We 
predict a trend of increasing planetesimal disk mass and 
decreasing compactness with stellar mass and metallicity 
and the lack of a giant planet in the system. However, 
these trends may be diluted due to the increased fre¬ 
quency of giant planets around more massive and metal- 
rich stars. On the whole, the range of planetesimal disk 
structures predicted by our model can account for the 
diversity of observed planetary systems. 
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4.2. Final assembly of planets 

The logical next step is to see what kind of planetary 
systems will assemble from the planetesimal disks formed 
in these simulations. This final step is, as of yet, an un¬ 
solved problem, but we will comment on the implications 
of our results. This work removes a large hurdle for the 
scenario in which STIPs form in situ. Our results show 
that high mass and compact planetesimal disks are not 
just possible but are the expected outcome of planetesi¬ 
mal formation and growth. Thus, if we ignore the effects 
of migration, STIPs will naturally form from in situ ac¬ 
cretion in these disks. 

If, however, we do not ignore migration, do STIPs still 
form? This is a more difficult question to answer be¬ 
cause simulations including migration have not relied on 
a massive inner disk and therefore did not include one. 
Additional work will need to be done to see if it is possible 
to grow planets farther out and migrate them through a 
massive planetesimal disk to form STIPs. Another pos¬ 
sibility is that planets accrete most of their mass as they 
move through the inner disk. 

Regardless of how the final assembly of planets occurs, 
an interesting question to be addressed in future work is: 
does the range of planetesimal disk surface density pro¬ 
files expected from these simulations produce planetary 
systems with the range of observed architectures? Ex¬ 
amining the more detailed statistics of systems produced 
from simulations may be our best bet in determining how 
planets are ultimately pieced together. 

5. SUMMARY AND CONCLUSIONS 

We have put forward a model for the growth of a plan¬ 
etesimal disk that is, in principle, straight-forward - plan¬ 
etesimal disks are built from the combination of direct 
planetesimal formation and the accretion of pebbles onto 
existing planetesimals. This model rests on two main 
assumptions: planetesimal formation occurs throughout 
the disk and pebbles rapidly grow in the outer disk and 
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then drift through the inner disk. 

In practice, our model is more complicated due to the 
complexity of pebble accretion and protoplanetary disk 
structure. This leads to a large number of model parame¬ 
ters, with large uncertainties, that affect the growth pat¬ 
tern of the planetesimal disk. However, one aspect that 
is common to our simulations regardless of these parame¬ 
ter values is that pebble accretion is efficient in the inner 
disk. This high efficiency allows for the formation of high 
mass inner disks without the need for a massive proto¬ 
planetary disk or enhanced planetesimal formation rates. 

The final structure of planetesimal disks in our sim¬ 
ulations is controlled by the rate of pebble delivery to 
the inner disk and the duration of this delivery. More 
massive protoplanetary disk and those with high dust to 
gas ratios will supply more pebbles to the inner disk and 
consquently build more massive planetesimal disks. The 
duration of pebble delivery is controlled by the lifetime 
of the protoplanetary disk and potentially the timing of 
giant planet formation. Giant planets will act as road 
blocks to any inward migrating pebbles thus starving the 
inner disk. With these comments in mind, we expect that 
more massive and/or metal-rich stars (with more massive 
disks and higher dust to gas ratios) and stars that do not 
host giant planets will contain more mass in terrestrial 
planets. 

Exactly how planetesimal disks evolve into planetary 
systems is still a matter of debate. However, our results 
strongly suggest that the initial conditions for late stage 
planet formation should include a massive planetesimal 
disk. This does support the idea that STIPs could form 
in situ without the need for the migration of material 
during or after planet assembly. However, this does not 
preclude the possibility that close in planets migrated 
significantly to their final orbits. 

The authors thank the referee for useful comments used 
to improve the paper. This material is based upon work 
supported by NASA under award No. NNX15AF02G. 


Hands, T. O., Alexander, R. D., & Dehnen, W. 2014, MNRAS, 
445, 749 

Hansen, B. M. S., & Murray, N. 2012, ApJ, 751, 158 
Johansen, A., Oishi, J. S., Low, M.-M. M., et al. 2007, Nature, 
448, 1022 

Johnson, J. A., Butler, R. P., Marcy, G. W., et al. 2007, ApJ, 
670, 833 

Kley, W., & Nelson, R. P. 2012, ARA&A, 50, 211 
Kretke, K. A., & Levison, H. F. 2014, AJ, 148, 109 
Lambrechts, M., & Johansen, A. 2012, A&A, 544, A32 
—. 2014, A&A, 572, A107 
Lithwick, Y., & Wu, Y. 2012, ApJL, 756, Lll 
Lynden-Bell, D., & Pringle, J. E. 1974, MNRAS, 168, 603 
Nakagawa, Y., Sekiya, M., & Hayashi, C. 1986, Icarus, 67, 375 
Ormel, C. W., & Cuzzi, J. N. 2007, A&A, 466, 413 
Ormel, C. W., & Klahr, H. H. 2010, A&A, 520, A43 
Pu, B., & Wu, Y. 2015, arXiv: 150205449 
Raymond, S. N., & Cossou, C. 2014, MNRASL, 440, Lll 
Santos, N. C., Israelian, G., & Mayor, M. 2001, A&A, 373, 1019 
Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337 
Terquem, C., & Papaloizou, J. C. B. 2007, ApJ, 654, 1110 
Testi, L., Birnstiel, T., Ricci, L., et al. 2014, in Protostars and 
Planets VI (University of Arizona Press), 339-361 
Valenti, J. A., & Fischer, D. A. 2005, ApJS, 159, 141 
Volk, K., & Gladman, B. 2015, arXiv: 150206558 
Wang, J., & Fischer, D. A. 2015, AJ, 149, 14 


10 


Weidenschilling, S. J. 1977, MNRAS, 180, 57 
Youdin, A. N., &; Lithwick, Y. 2007, Icarus, 192, 588 


